Glass Transitions and Shear Thickening Suspension Rheology 

C. B. Holmes and M. E. Catefl 
School of Physics, The University of Edinburgh, JCMB, 

S : "' ''""'' ''"""'"'"• "''"'""'"'• ^"^ ""' ""■ 

^ ■ M. Fuchs 

C ■ 

^ ' Fachbereich Physik, Universitdt Konstanz, D- 78457 Konstanz, Germany. 

>—i ; 

00 ; p. SoUich 

I— I, Department of Mathematics, King's College, 

^^ \ University of London, Strand, London, WC2R 2LS, UK. 

c/3 ; (Dated: February 2, 2008) 

C^ . We introduce a class of simple models for shear thickening and/ or 'jamming' in 

^ \ colloidal suspensions. These are based on schematic mode coupling theory (MCT) 

i-rt I of the glass transition, having a memory term that depends on a density variable, 

C ' and on both the shear stress and the shear rate. (Tensorial aspects of the rheol- 

Q ■ ogy, such as normal stresses, are ignored for simplicity.) We calculate steady-state 

'""'. flow curves and correlation functions. Depending on model parameters, we find a 

T-H ! range of rheological behaviours, including 'S-shaped' flow curves, indicating discon- 

J^ \ tinuous shear thickening, and stress-induced transitions from a fluid to a nonergodic 

^S| ' (jammed) state, showing zero flow rate in an interval of applied stress. The shear 

''^ ■ thickening and jamming scenarios that we explore appear broadly consistent with 

^^ . experiments on dense colloids close to the glass transition, despite the fact that we 

■^ , ignore hydrodynamic interactions. In particular, the jamming transition we pro- 

--^ \ pose is conceptually quite different from various hydrodynamic mechanisms of shear 

d ' thickening in the literature, although the latter might remain pertinent at lower 

Ch ' colloid densities. Our jammed state is a stress-induced glass, but its nonergodicity 

i-rt ' transitions have an analytical structure distinct from that of the conventional MCT 

Ch , glass transition. 
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INTRODUCTION 

Under flow, concentrated hard sphere suspensi ons can shear thicken : the measured viscos- 



i ty inc r eases with the app l ied stress or flow rate iBender and Waenerl ( 19961 ) , Id'Haene et al. 
( 1993[ ). lFrith et al\ (jl996t ). lLaunl (jl994^ ■ l0'Brien and Mackavl (J2nnnl )]. The effect can be pro- 



nounced, in some cases giving an order-of-magnitude jump in the shear stress on exceeding a 
'critical' shear rate. In such 'discontinuous' shear thickening the flow often become s erratic 



sugges ting the continual formation and breakup of stress-supporting structures (Frith et al. 



( 19961 )]. This idea is supported by recent pressure-driven flow experiments on model colloids. 
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in which transient sohd regions are directly visuahsed Hawi (200J)]. Computer simulations 
give a similar picture, with simulations at infinite Peclet numb er suggesting that , in th e 
absence of Brownian motion, a static, load bearing jam may form Melrose and Balll ( 19951 )]. 
This is reminiscient of jamming in gran ular materials, in w hich the solidity of the material 



may result from the applied stress itself Gates et a,l\ ()l998t )]. 



The notion that strongly shear thickened st ates have the characte r of an amorphous solid 



or glass, is supported by recent experiments jBertrand et al\ (J2| 



in which a droplet of 

shear thickening material undergoes a fiuid-solid transition upon shearing. In these exper- 
iments, a suspension (mean particle radius a ^ 3/im) was transformed into a solid state, 
simply by stirring with a spatula. The jammed solid persists after cessation of stirring, but 
is refiuidised by vibration. Also, it is reported that introducing a droplet of the quiescent 
liquid state adjacent to a jammed droplet (at the same concentration) results in the two 
coalescing to form a l arger fluid droplet. The description of these phenomena reported in 
Bertrand et al\ ( 20021 ) is brief, but similar observations have been made in suspensions of 
zeolite particles and hard-sphere colloids . (Note that the system of lBertrand et all ( 20021 ) is 
not hard-sphere like, and the reported phenomena are at relatively modest concentrations. 
For hard spheres, the relevant concentration range lies above 50% by volume.) Thus there 
is some evidence of two distinct states of a dense suspension, one fluid and one jammed. 
These must presumably have different internal organization. Indeed, the surface of jammed 
droplets appears dry or dull (becoming shiny again only when refiuidised by vibration or 
coalescence) suggesting that partial emergence of particles is needed to sustain the jammed 
state. This is consistent with the presence of large, static deviatoric stresses within the bulk 
of a jammed droplet, balanced at the surface by capillary forces. 

Such experiments pose important challenges to the rheological modelling of colloidal sus- 
pensions. In particular they suggest, though they do not prove, that static shear stresses 
might be sustained indefinitely in a bulk sample at rest in a rheometer - even if the 
same sample is fiuid at lower stresses. (In this paper we address only this bulk geome- 
try; we hope to return to freestanding droplets in future work.) On the other hand, shea r 
thickening is commonly a ttributed to ' h ydrodynamic clustering' |Ball and Melrosd ( 19951 ). 
Bergenholtz et al\ ()2002l) . JBradvl ()l99fi[) . JMehose and Balll ()l995[) ] - that is, the formation 
of transient clusters of particles in lubrication contact. This results in large hydrodynamic 
stresses in fast fiowing suspensions, and there is simula tion and other evidence th a t these 
are th e dominant str e sses in the shear-thickened regime Bender and Wagnerl ( 19961 ). iBradv 
( 19961 ). IPhung et al\ ( 19961 )]. However, some of these simulations are limited to relatively 
modest volume fractions of colloid (0 < 0.49); experiments at much highe r concentrations 
suggest instead a thermodynamic mechanism |0'Brien and Mackav (120001)] . In a ny case it 
seems highly unlikely that the 'static jamming' observed bv lBertrand et all ( 20021 ) is due to 
hydrodynamic interactions alone, since in the persisting solid phase there is no macroscopic 
fiow to provide those interactions. Also, as far as we know, hydrodynamic models for shear 
thickening offer no immediate exp lanation of nonmonotonic regions o f the fl ow curve. Such 
regions may be observed directly Bertrand et a/.l ( 20021). Frith et all ( 19961 )]. and can lead 
(via well understood shear-banding mec hanism s Olmstedl (|l999r i]) to discontinuous shear 
thickening and accompanying hysteresis LaunI (199^]. 

Thus it seems clear that mechanisms othe r than pure hydrodynamics are at work in 
sufficiently dense shear-t hickening suspensions lO'Brien and Mackav^ (120001 )]. Indeed, it has 

u 



already b^en empha.i.ed | Rall and M.>osj j lPJtA Ulros. a.d Ra lU SH) that deviation. 



from pure lubrication forces can dominate the physics of any hydrodynamically clustered 



state. (These deviations could arise, for examp l e, from deform a tion o f stabilizing poly- 
mer layers around the parti cles Ball and Melrosd ( 19951 1 , iMelrosd ( 20031 ) , iMelrose and Balll 
( 19951 ) . iMelrose et al\ ( 19961 )] . ) Cruc i ally, h ydrodynamic theories neglect something else too: 
the colloidal glass transition Puseyl (|l989)]. 

In quiescent colloidal fluids, Brownian relaxation processes are often interpreted in terms 
of particle 'caging'. There is a short-time relaxation, called /3-relaxation, which is a many- 
body process during which particles remain within a cage of nearest neighbours. This is 
followed by a late-time a-relaxation, which corresponds to particles escaping their cages. 
As the glass transition is approached, the cages become tighter, and both of these relaxation 
times become longer. In quiescent dense suspensions, it is this mechanism (rather than 
simply the vanishing of mobility as interparticle gaps close to zero separation at random 
close packing) that leads to the observed slow relaxations and large viscosities. Within 
the so-called mode coupling theory (MCT) of the glass transition, the divergence of the 
a-relaxation time (and hence the viscosity) occurs at a colloid volume fraction 0^ around 
52% (although this value is very sensitive to the input static structure factor). This is 
far below the random close packing density (64%) which is the natural location of any 
divergence in purely hydrodynamic theories. In practice, this MCT value is too low, but 
not by much - light scattering d ata on a-relaxation locate the glass transition at 58% 



van Meeen and UnderwoodI ()199 



We contend that the slow relaxation of thermodynamic forces caused by particle caging 
remains important in sheared dense suspensions, and is key to understanding their non- 
Newtonian rheology. If this is true, the interesting physics should arise at shear rates 7 set 
by a renormalized Peclet number 7X0 involving the a-relaxation time. In contrast, hydrody- 
namic mechanisms are controlled, in the simplest approximatior i, by the bare Peclet numbe r 
7rn, with To the diffusive relaxation time of an isolated particle Bender and Wagneii ( 19961 ). 
Phung et ai, (1996) ]. Less naively, one could expect this to be renor malized to 7T1 where T ^ 
is a lubrication-diffus i on tirn e that diverges at random close packing Foss and Bradvl ( 20001 ) , 
O'Brien and Mackavl (|2000( )]. But as emphasised above, even ri is vastly shorter than Tq, in 
fluids close to the glass transition, and infinitely shorter throughout the glass phase itself, 
where r^ has diverged. 

While in practice the hydrodynamic and thermodynamic forces in very dense colloids 
may be interdependent, this separation of time scales makes it legitimate to explore the pos- 
sible connection between the glass transition and shear thickening in its own right, ignoring 
hydrodynamic interac tions in the first inst ance. From this viewpoint, the jammed state 
found experimentally pertrand et al\ (120021 )] is a candidate for a stress-induced, anisotropic 
colloidal glass. Such a state could arise if the applied stress alters the material's struc- 
ture, leading to an increased number of close contacts, hindering diffusion and promoting 
dynamical arrest through the tightening of the cage around each particle. 

In this picture, as in a conventional glass, solidity is due to thermodynamic forces (purely 
entropic for hard spheres): deformation of the material results in a free energy penalty which 
appears macroscopically as elasticity. An alternative is that Brownian motion is sufficiently 
weak that these entropic forces are irrelevant, leaving interparticle forces alone to resist 
deformation: in hard spheres, stress applied to the jammed solid is then transmitted by 
direct contacts. In this case, the jamming more closely r esembles that of a gra nular medium 
|Cates et all ( 19981 )]. (In fact it has been suggested bv iLiu and Nagell ( 19981 ) that, even in 
this non-Brownian case, arrest is related to a glass transition; we do not pursue this here.) In 
a jammed suspension, it is not clear which mechanism dominates; the answer might depend 



strongly on particle size. Nor is it clear under what conditions static jamming may be 
achieved: sufficient Brownian motion might give all or most iammed state s a finite lifetime. 
Note that, in contrast to some authors (for example iLiu and Nagell ( 199a l). in this work we 
use the word 'jamming' to mean an arrest transition caused specifically by imposed stress 
rather than other causes, such as attractive interactions. 

Our modelling strategy, which was outlined with selected results in lHolmes et al\ ()2nn3r ). 
neglects hydrodynamics, and assumes a dynamical gl ass transition pictu re rather than one 
based on direct mechanical contact between particles Gates et al\ (|l998t )]. It builds upon a 
schematic version of MCT, and is guided by more formal ( and less sch ematic) MCT work 
on the role of shear in the glass transition Fuchs and CatesI (|2002ll2003h ]. The latter formal 
approach currently gives only shear-thinning behaviour, or strain-induced fiuidisation. To 
this we add a new feature in our schematic models: stress-induced arrest. (Capturing the 
same physics, beyond t he sche r natic level, remains an open issue for future development of 
the more formal theory CatesI ( 20031 )].) In combination, these two features suggest various 
shear thickening and jamming scenarios, some of which appear broadly consistent with the 
experimental features noted above. 

In some but not all cases, these scenarios are similar to ones found before from a quite 
different approach to glassy rheology, originally intend ed for 'soft glasses' (foams, dense emul- 
sions) rather than colloidal suspensions. In particular. iHead et al\ (|200lL 120021) incorporated 
an ad- hoc local jamm ing effec t into the 'soft glassy rheology model' (SGR) Cates and SoUichl 
( 20041 ) . iFielding et al. (20001 ). ISoUich et al\ ( 1997 )] thereby converting its shear-thinning be- 
haviour into a shear thickening one. Our work on schematic MCT, though equally ad-hoc 
at this stage, has the advantage that one day it might be underpinned by full MCT-type 
calculations, whose form might well be guided by comparing our different schematic model 
variants with experimental data. (A drawback is that, unlike SGR-based work, we can only 
address steady state properties a nd not dynamical one s.) Some of the flow curves we obtain 
are also related to those found bv lHess and Hesa ( 19941 ) although their work did not attempt 
to relate shear thickening to glass transitions. 

The rest of this paper is organized as follows. Section ^1 provides a brief introduction to 
MCT; in Sec. 11111 we formulate our new class of schematic models. Analytical and numerical 
methods are described in Sec. HVland the results of these presented in Sees. El and EH Sec. 
IVllI contains a discussion of our results in the context of experiment, and we conclude in 
Sec. EnH 



II. MODE COUPLING THEORY 



For a comprehensive description of MCT, see iGotzd ( 19891 ). Here, we outline only those 
aspects pertinent to the current work. The central quantities of the theory are the Fourier- 
space density fluctuations at wavevector q, 6p{q,t). The correlators of these qua ntities, 
(i)r ^(t) = (5p(q, t)5p(— q, 0))/(|(5p(q)P) may be measured in scattering experiments PusevI 
()l989l )]. and describe the system's dynamics. In a liquid, the system is ergodic and 0q(t) 
decays to zero with time for all q. In a glass it does not: limt^oo 0q(^) = /q > 0, where the 
nonergodicity parameters /q characterise the arrest. A flnite /q implies inability to relax on 
a lengthscale ~ 27r/g. (In general, all /q become nonzero at the same transition point; once 
one Fourier component of density is frozen, the others see this as a random potential and 
acquire non-decaying mean values.) 

Within MCT, equations of motion are found for the correlators 0q(t). Assuming over- 



damp ed local motion appropriate to colloids, and dropping q-subscripts, the result is Gotzd 
(19981)]: 



0(t) + To^{t) + [ m{t- t')^{t') dt' = 0, (1) 

Jo 

where To is a timescale set by the microscopic dynamics. Here m{t — t') is the memory 
function, and describes a retarded friction which, in the colloidal glass transition, arises by 
caging of a particle by its neighbours. In MCT, the memory function is found approximately 
by integrating (over wavevectors) a quadratic product of correlators, with coupling constants 
that depend solely on the static structure factor S{q) = (l/iV)(|5p(q)p) of the system. The 
To term represents instantaneous solvent friction (with many body hydrodynamics omitted) 
and gives exponential relaxation of correlations in dilute colloids, in which memory effects 
are negligible. Note that the nature of the interparticle interactions enters through S{q); 
MCT is not limited to hard-sphere colloids. 

In MCT, the zero shear viscosity diverges with the a-relaxation time of the density corre- 
lations. An infinitesimal stress applied to a colloidal glass state gives zero steady-state shear 
rate 7: glasses are solid, although creep (sublinear in time) is not excluded. The extension 
of MCT to finite 7 requires care, as Eq. [T] relies upon the fiuctuati on dissipation theo- 



rem (F DT) jBouchaud et all (|1996[ )] which need not hold under shear |Barrat and Berthier 



( 20001 )]. However, extended forms of MCT d escribing shear have been formu lated, either 



wi th the assumption of an FDT- like condition JMivazaki and ReichmannI ( 2002i )]. or without 



it 



Fuchs and Cated (|2002ll2003t )]. These approaches lead to equations of motion for correla- 



tors like 0q (redefined to account for simple advection) which formally resemble those of the 
unsheared theory. Solution of these equations, performed via MCT-like approximations, 
leads to the conclusion that structural relaxation always speeds up in a fiowing system: 
Tni'i) < Tn(0). (Ana logous behaviour is seen in mean field spin models of driven glasses 
Berthier et al\ ( 20001 )].) Physically, this is because small scale diffusive motions are magni- 



fied by advection, facilitating break up of the cage on a timescale set by I/7. In q-space, 
the important fiuctuations, which have wavevectors near the peak of S{q), are advected on 
this timescale to higher q where t hey decay rapidly . Even those fluctuations which are not 
directly advected become ergodic Fuchs and Cated ( 20021120031 )]. 



So far then, these microscopic extensions of MCT to sheared systems do not predict jam- 
ming under shear stress. This may be because shear thickening requires hard-core particle 
interactions which are not fully captured by the harmonic approximation to the thermo- 
dynamic forces conventionally made in MCT Cated ( 2003^ ]. The latter is tantamount to 



an expansion of the free energy to second order in density fluctuations; in order to see the 
jamming effect, one might have to go beyond this order. We do not attempt this formidable 
task here: rather, we adopt a simpler approach, beginning from 'schematic models' of the 
MCT glass transition. These stripped-down models represent a gross simpliflcation of the 
full t heory and ye t manage to capture many of the key features of the quiescent glass transi- 
tion Gotzd (|l989l )]. In them, one considers a single correlator 0(t) (which may be considered 



to represent fluctuations at some typical wavevector), rather than the inflnite set 0q(t). The 
memory function is then written as a polynomial of this correlator, with coefficients (cou- 
pling constants) that encode the interactions between particles (just as S{q) encodes these 
in full MCT). Larger coupling constants correspond to higher colloid densities and/or lower 
temperatures. 

Of interest here are two schematic models, known as the F2- and F12-models, defined 



6 

by a shared equation (Eq. 2]) for the correlation function 0(t), and by the memory functions 

m{t) = f20^(t) (F2 model) (2) 

m(t) = fi0(t)+t;20^(t) (F12 model) (3) 

Each model exhibits a glass transition with increasing f„; fo r the F 2 model, this occurs 
at f 2 = 4 whilst in the F12 model (which has f 2 > 1 |Gotzd (19891)]), there is a locus of 



transition points given by "Wi = f 2 ( A/4/f2 — 1 ) • As these transitions are crossed, the single 

nonergodicity parameter / = 0(oo) jumps discontinuously, from zero on the liquid side of 
the transition, to a finite value fc on the glass side. On further increasing the coupling, 
there is a nonanalytic (square root) increase in / beyond fc- This behaviour of / is faithful 
to the full MCT in which a similar discontinuity and square root law are found for the /q. 

III. RHEOLOGICAL MODELS 

To incorporate shearing into the schematic MCT models we need to address two distinct 
effects. Firstly flow erodes memory: particles with nonzero separation in the flow gradient 
direction become separated in the flow direction on times t ~ I/7. Thus, diffusive motion 
is accentuated by the flow, with the result that cages are broken and correlations decay 
to zero. The second effect is that stress (as distinct from flow) can alter the local caging 
of particles such as to promote arrest. Within a schematic theory, the latter requires us 
to postulate a stress dependence of the coupling constants. Perhaps the simplest idea is 
to assume a Taylor expansion of the coupling constants Vn {n = 1, 2) in Eqs. \2fS\ that is, 
Vn{cr) = Vn + J2k'^kO'^- ^^^ could then neglect odd powers on the grounds of symmetry. 
However, this approach is limited to the region around cr = 0, whereas we are interested in 
nonzero stresses (in particular, values close to the yield stress of an incipient glass). Also, 
it turns out that the topology of flow curves depends strongly on the form of Vn{cr) at large 
a, rather than that near the origin. Hence, although we choose mainly simple power laws 
for f„(o"), these are not to be viewed as a Taylor expansion around a = 0. In particular, a 
linear dependence w„(cr) = f „ + aa, is not excluded by symmetry. Here, and from now on, 
a denotes the absolute value of the shear stress. The various model forms chosen for Vn{cr) 
are specified below in Eqs. ITHTTl 

To allow for the flow-induc ed memory loss, the memory fu nction should also be attenuated 
by a strain dependent factor Fuchs and Gated (|2002ll2003l ). This suggests that we write 



m{t)=Y,Vn{a)r{t)fm, (4) 



where the fu nction f{'jt) is monotonically decreasing. Microscopic theory Fuchs and Gated 



( 2OO2L 120031 ) suggests a specific form: f{'yt) = (1 + 7^t^) ^. But for algebraic convenience 



we shall mainly consider an exponential form fXjt) = exp (— 7^). 

In Eq. 0] for the memory function, the effect of shearing appears in two quite separate 
ways: via the applied stress a, and with the resulting (steady state) shear rate 7. These two 
quantities are related through the viscosity 

V = ^/7- (5) 

In order to close the model, we now need a prescription for the viscosity. In linear re- 
sponse theory, the shear viscosity is expressed as the time integral of a stress correlator 



Hansen and McDonaldl ( 19761 )]. In our simplified model, we have but one correlator, and 



so we choose schematically to identify the viscosity as 

r]= / 0(t)dt = r, (6) 



where the appropriate elastic modulus has been set equal to unity, rendering the stress 

dimensionless. Eq. IHl is somewhat ad hoc, since we are not studying the linear response 

regime, but experim ent suggests a close relationship between viscosity and relaxation times 

Meeker et al\ (1997|)]. Eq. IHl means that dynamical arrest (nonzero /) implies a divergent 



viscosity, and hence zero shear rate, in steady state. This is fully consistent with arguments 
made above that finite 7 precludes arrest. Note that equating t] to an integral of 0^(t) (say) 
would give very similar results. 

Eqs. ^E]and|niare shared by all the models in the class we study here. Different members 
of the class have different choices of Vn{cr) and f{'jt) in Eq. HJ Among the variants we shall 
address are 

m{t) = (f2 + a(T)02(t) exp (-7^) Model I (7) 

m{t) = {v2 + aa^)0^(t) exp (-7^) Model II (8) 

These variants are amenable to an analytic approach, which we present in the next Section. 
We also consider a Model III which interpolates between Model II at small stresses and 
Model I at large ones; this has 

m{t) = {v2 + ainf [a, a^])(j)'^{t) exp (-7^) Model III (9) 

The other variations which we shall study are: 

m{t) = {v2 + aa)02(t)/(l + ft'^) Model IV (10) 

m{t) = [{vi + aia)(f){t) + {v2 + a2fT)0^(t)] exp (-7t) Model V (11) 

Model IV differs from mod el I only in that the choice o f f{'jt) is closer to the form ex- 
pected from the full MCT JFuchs and CatesI 1)20021 12003[ )]. Unfortunately this limits ana- 



lytic progress, but we present below numerical results for this case. Model V departs from 
variants I-IV by being based on the F12 model rather than F2; similar remarks apply. 

In each of these model variants, there are control parameters such as f„ and a: these 
are expected to have nontrivial dependence on both the concentration of colloids and their 
interactions. Since in practice a range of different shear-thickening scenarios are observed for 
different materials, an attempt at a parameter-free theory would, at this stage, be misguided. 
We choose o; > in this paper. 

IV. SOLVING THE MODELS 

A. Analytical Developments 

Our analytical work is based on the F2 model and hence restricted to Models I-IV, for 
which vi = 0. The methods are adapted from standard techniques used in solution of the 



MCT equations of motion (see, eg Gotzd (19891)]). We begin by defining a Laplace transform 



f{z) = i / f{t) exp{izt) dt; Q{z) > 0. (12) 

Jo 

Then, Laplace transforming Eq. ^ we obtain 

-1 



z- 1 



(13) 



iTo+m{z) 



The dynamical slowing down associated with the glass transition is contained in m{z). In 
contrast, the instantaneous friction term remains constant, and so, close to the transition, 
\to\ <^ |'^(-2)| for low frequencies. For clarity in the following, we label the position in 
parameter space by the vector v = (f 2, a, cr) so that, eg, (l)v{z) is the frequency-space corre- 
lator at frequency z for parameters given by v. (Likewise, 7v is the shear rate for the given 
parameter values.) At low frequencies then, we neglect the instantaneous friction term to 
give 

^'^-^^^ = ,m (z) (U) 

We now split the correlator into the sum of its late time limit and a (time-dependent) 
remainder, writing 4>v{t) = fv+gv(t). Now we assume that, for suitable values of the relevant 
parameters, our model admits nonergodic solutions. In a nonergodic state the relaxation 
time T is divergent and so 7 = for any finite stress. Therefore the zero-frequency limit of 
Eq. HH gives 

Vfl (15) 



l-/v 

where V = ^2(0"). Here, we have used the fact that, since gv(t) tends to zero at late times, 
terms o f order Qv(t) or q^ft ) give rise to singularities weaker than I/2; as 2 — *> 0. It has been 
shown Gotze and SjogrenI ()l995r )] that the nonergodicity parameter (limt^oo 0v(^)) is the 



largest of the solutions to Eq. IT^ 

This result fixes the nonergodicity parameter of an arrested state with coupling V: there 
are nonzero solutions for /v provided V^ > 4, so that a glass becomes possible beyond this 
value of the coupling. However, this does not imply that we necessarily have a glass transition 
at V^ = 4, since we have calculated this by assuming arrest, which requires 7 = 0. (Without 
shear, our model is simply an F2 model for which the result V"c = 4 is well known.) In Sec. 
im we stated that passing through a static MCT glass transition entails a discontinuity in 
the nonergodicity parameter / followed by a nonanalytic variation as the coupling constants 
increase. Solving Eq. ^J shows that such behaviour will arise in Models I- IV only for an 
arrest transition occurring at the F2 value of V" = 4. We shall return to this point later. 

1. Asymptotics for correlator and memory function 

We now consider the behaviour of correlation functions close to a transition at (say) 
Vc, corresponding to some V > 4. Near such a transition, |v — v^ is a small quantity. 
To identify other quantities which may be treated perturbatively, we recall static MCT 
behaviour: nonergodic solutions to Eq. [T] show a decay, at the /5-timescale r/3, onto a 
plateau, with 0v(t > r^) ~ /v For V > 4, on timescales for which %t <^ 1, the difference 



4>v{t) — (f)v^{t) is also a small quantity. Provided 7v vanishes smoothly at the jamming 
transition, we can ensure that this condition holds at and beyond the timescale rg on which 
4>v^(t) reaches its long time limit. That is, we choose |v — vd such that JTfs -C 1. (Our 
ability to do this depends on the fact that rg is finite, unlike t^, whenever V > 4.) Thus we 
have identified a time regime where gv{t) = 4'v{t) — f^^ is a small quantity. In the fre quency 
domain, this gives z(j)v{z) = —fv^ + zg^iz). In line with standard MCT techniques [Gotzel 
f 19891 )] ■ we now treat gv{z) perturbatively in the Laplace transformed equation of motion. 



Eq. HH Expanding the LHS in powers of zg^{z), we find 

z(t)w{,z) -/v, , zg^{z) z'^gKz 



1....W i-;,/ir^^-ir"t].-''(^^«vW). (16) 

We now wish to approximate the RHS of Eq. Elin the same regime of frequency. By consid- 
ering the memory function in the time domain, j^t and gv{t) can be treated perturbatively. 
So can the distance from the transition in parameter space: we write V = Vc + e where Vc 
is the value oi V = ^2(0") at Vc. 

To make further analytic progress we now restrict the form of the function f{'yt) in Eq|3] 
to be exponential. (This excludes Model IV; we have already excluded Model V.) We find 
the memory function then to obey 

mv(t) = (K + e) {fl + 2UMt) + 9l{t)) (l - 7vt + 0(7v^')) • (17) 

In the time regime of interest, terms 0{gv(t)%t) and O {gl{t)) are negligible with respect to 
terms O {gy,{t)). Additionally, terms O {eg^{t)) and O (e7vt) may be neglected with respect 
to terms 0{t). The RHS of Eq. IT^then becomes 

zm^{z) ^ -Vfl + 2Uy,zg^{z) + VJlj^t/z. (18) 

Using Eq. ^J we find that close to the transition Eq. IT^ becomes 

zg.X^) (^r^T^y - 2/vc^c) - VJl^.i/z + fie = 0. (19) 

Inverting the Laplace transform, we find the correlator in the time domain to be 

or, expressed more succinctly, 

0v(t) = /(l-t/f) (21) 

with nonergodicity parameter / and a relaxation time f defined appropriately. 

Eq. 1^ is consistent with an exponential decay from the plateau at late times, (f) ~ 
/exp[— t/f], a result also found numerically in Sec. IVl below. In what follows we further 
assume that this is the correct late-time form for the correlator; this allows us to find 
explicitly the locus of arrest transitions arising in Models I-III. 

We begin by evaluating Eq. Uni at 2; = to deduce the existence of a terminal relaxation 
time 

T = To- im{0). (22) 
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As we approach a jamming transition at finite stress, the relaxation time f defined in Eq. 
ITT] must diverge (since the shear rate vanishes at the transition point). We continue to 
assume that, at Vc, V^ > 4 in which case the r^ remains finite. Therefore, approaching the 
transition, the terminal time r is primarily set by the late-time exponential relaxation of 
(pv(t). Accordingly, Eq. 0^ becomes 

T = To + V exp(-7t)02(t)dt (23) 







POO 

kvI exp(-7t)/^exp(-2t/f)dt, (24) 

Jo 

where we have approximated the correlator by 0(t) = /exp(— t/f) over the whole temporal 
range (neglecting deviations at t < rg), and dropped a (non-divergent) contribution from 
the regular part of the memory function. By definition, r = J°° 0(t) dt , and so we require 
/o°° /exp(-t/f) = r, giving f = r//. 

We now write Eq. |^ as an iteration: the nth approximation to the relaxation time r is 
calculated using the (n — l)th value to determine the shear rate. Thus 



r' 



/•oo 

-(^) ^ / V/^e-'^^/^We-^/V-CDdt, (25) 

Jo 

which we can solve to give the ratio of r^^^^ to r'^^^: 

:^ = ^{Vf-2). (26) 

Crucially, this ratio is independent of r^^^: thus, no matter how large the initial guess, the 
iteration leads to a larger value if the RHS is greater than unity. Where this is the case, the 
relaxation time is divergent. Thus, setting the RHS in Eq. |2nito unity locates the arrest 
transition. Recalling the definition V = f2(cr), this lies at a stress Uc which obeys 

fc[v2{cTc)fc-2]-a, = 0, (27) 

where fc is given by the largest solution of /c/(l — fc) = V2{(Jc)fl (from Eq. IT^ . 

Except for the special case ac = 0, these transition points do not coincide with the F2 
transition at \^ = 4, and have a different character to MCT transitions in general. From a 
rheological point of view, the transitions give values of the stress for which the shear rate 
first becomes zero: these are jamming transitions. Their locations according to Eq. 1771 will 
later be checked against the numerical solutions for the full fiow curves. 

2. Large stress limit 

We turn now to investigate the behaviour in the limit of large stress. We begin by 
considering cases in which the relaxation time diverges in that limit: in this case, Eq. QUI 
implies 

lim afv2{(r)/a > 1. (28) 

cr— >00 

The implications of this result depend upon the model: if f2(cr) increases sublinearly at 
large a, this criterion will not be satisfied for any finite a and so fiuid-like behaviour is 
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always recovered in the limit of large stress. For super-linear V2{cr), as arises in Model II, 
this condition places no restriction upon a beyond a > 0. Hence Model II always has zero 
shear rate at high stresses: the jammed state remains solid indefinitely. 

The most interesting case is when ^2(0") increases linearly at large a, which applies in 
both Model I and Model III. Here Eq. EHlcan only be satisfied if a > 1, since in a jammed 
state the nonergodicity parameter / will approach unity in the limit of large stress (see Eq. 
ITH|) . This implies that, for a < 1, the shear rate can be made arbitrarily large by increasing 
the stress. In this case, we can calculate the relaxation time in the large stress limit from 
Eq. |2ni Only the region t ~ I/7 contributes significantly to the integral, and for large shear 
rate s this r e gion becomes arbitrarily small. Thus, we may expand the correlator in powers 
of t JHinchI fll99l| )]. giving 



/ e-'^'^'^{t)dt^ / e-^*(02(O) + O(t))dt. (29) 

Jo Jo 

The initial condition 0(0) = 1 (normalisation) gives 

/"OO 

T^To + V2ia) / e-'^^/^dt (30) 

Jo 

which is solved in the limit of large stress by 

r=^ (a.<l). (31) 

If r grows sublinearly with a at large stress, both quantities diverge in this limit. In this case 
a > 1 and r = ^^i demanding a non- negative relaxation time leads to the conclusion that 
a = 1. Thus, in Model I, a = 1 is a watershed between ergodic and nonergodic solutions in 
the limit of large stress. For a < 1, solutions are ergodic with the viscosity approaching a 
limiting value; but for a > 1, the flow rate 7 is zero at high stress. 

B. Numerical Solutions 

The preceding analytical results are useful but restricted to certain model variants. Also 
they do not provide closed form expressions for flow curves (only transition points). We 
have therefore also solved ou r schematic rheolog ical models numerically, using an established 



algorithm for MCT models [Fuchs et ali (|1991|)]. The numerical method involves re- writing 



the governing equation in such a way that 0(t) is determined from its value at earlier 
times, allowing iterative calculation of correlation functions. Numerical efficiency is aided 
by treating slowly varying quantities as constant over the timestep. To adapt this algorithm 
to our models, we employ the following procedure: 

1. Input the parameters of interest. 

2. Make an initial guess at the relaxation time r for these parameters. 

3. Calculate the shear rate using this relaxation time and the chosen value of the stress. 
This defines iteration-dependent coupling constants fn(o") and the memory loss func- 
tion /(7t). 
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4. Solve the resulting model via the methods of lFuchs et all ( 19911 ) (holding the coupling 
constants fixed). 

5. Re-calculate the relaxation time by integrating the resulting correlator. 

6. Repeat from step 3 onwards using the new value of the relaxation time until converged. 
The result is a self-consistent solution. 

Numerical solutions using this algorithm were obtained using a FORTRAN code. We 
checked that the code returned the expected F2, F12 results for a = 0. For some parameters 
close to a nonergodicity transition, the code became slow to converge, especially for 'bad' 
initial guesses. This problem was circumvented by starting the iteration off with a range of 
initial values and bracketing the correct value according to whether the initial trajectories 
showed increasing or decreasing r. In some cases this generated improved initial values 
from which the above algorithm could attain convergence; in others, the process was simply 
continued until converged to the required accuracy in its own right. (The stability of results 
obtained in this manner was then checked by a final iteration.) 



V. RESULTS: MODEL I 

In this Section we present detailed results for Model I. Later fSection lVI|) we shall examine 
the effect of altering the form of the memory function, resulting in Models II-V. (This is 
easier than attempting to present results for all the model variants in parallel; it does not 
necessarily mean that we prefer Model I to the rest.) From now on we use units where 
the local relaxation time far from the glass transition, Tq, is unity. We continue to set the 
transient elastic modulus G to unity and can therefore identify the viscosity rj with the 
terminal relaxation time r. 



A. Steady-State Rheology 

1. Choice of parameters 

We expect the most interesting behaviour to arise close to the quiescent stress-free glass 
transition, where there is maximal interplay between the long relaxation times and the 
effect of shearing. For hard-sphere colloids, discontinuous thickening tends to occur for 
concent rations (/> > . 5, for which the zero-shear viscosity is around 50 times that of the 
solvent JMeeker et all ( 19971 )]. Thus, we first study values of V2 for which t > 50 (recall 



that To = 1). The stress-coupling parameter a does not have such an obvious experimental 
analogue, and it is less clear what values to investigate. However, from Sec. IIVI we expect 
values of order unity to lead to stress-induced nonergodicity, so we mainly consider such 
values. (As indicated in Sec. IIV A21 and confirmed below, a = 1 in fact separates differing 
regimes of behaviour.) Stress can result in thickening, in principle, whenever V2 + acr ~ 4: 
but whether this actually occurs depends on the interplay with shearing- induced memory 
loss. 
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FIG. 1: Model I flow curves for a = 0.95: points are numerically determined values of the shear 
rate for a given stress. Shear thinning is seen and, for vi = 3.8, thickening is also discernable at 
low stresses. 



2. Dependence on V2 

We begin by varying V2 at a fixed value of a, chosen slightly less than unity. Fig. ^ shows 
the flow curves for a = 0.95 with 3.5 < f 2 < 3.8. In all cases, there is a regime of shear 
thinning and, for the largest value of V2, the shear thinning is preceeded by thickening at 
small applied stresses. In all flow curves presented (here and below) lines are to guide the 
eye. For more modest values of V2 the rheology, not shown, is close to Newtonian. (For 
V2 = 3, for example, the viscosity changes only by about 15% upon increasing the stress 
between zero and a = 5.) 

Fig. 121 shows how the trend towards shear thickening develops upon further increasing 
V2- In contrast to the rather gentle curves of Fig. ^ the behaviour observed here is rather 
drastic. As V2 is increased above 3.8, the shear thinning behaviour shown in Fig. ^ is 
interrupted by shear thickening whose severity ranges from modest to extreme - for the 
larger values of f 2 shown, the thickening becomes so strong that regions of negative slope 
appear. In these regions, increasing the stress lowers the shear rate. S u ch regions of neg ative 



appear, in these regions, increasing the stress Lowers the shear rate, a u ch regions 01 neg ative 
slope should be mechanically unstable to shear banding DhontI ( 19991 ) . lOlmstedl (|l999l )] and 
should lead in practice to discontinuous shear thickening as we discuss in Sec. IVllI 

For the largest values of V2 shown, there is a window of stress for which the viscosity has 
diverg ed, resulting in "full j a mming" - the creatio n of a nonflowing state by application of 
stress JHess and Hesd ( 1994 ). iHolmes et all ( 20031 )]. Within this window, the only solution 
is a nonergodic, jammed state. The edges of this window are the jamming transitions whose 
positions are given analytically by Eq. |2ZI The analytic and numerical results for their 
positions are in agreement. We examine in more detail the nature of these transitions in 

Sec. En 

The above scenario describes the stress induced arrest of a state with V2 < 4; this is a 
liquid at low stresses and therefore the region of the flow curve close to the origin is always 
Newtonian (albeit with a viscosity that diverges as f2 — > 4). By setting V2 > 4, we can 
consider shearing a system that is arrested in the quiescent state (the usual definition of a 
colloidal glass). In this case, as shown in Fig. El there is no steady-state flow unless the 
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FIG. 2: Further Model I flow curves for a = 0.95. For the two largest values of V2, there is a 
window of stress for which the viscosity diverges. The limits of this window (denoted ad and ac2, 
as shown here for V2 = 3.87), are jamming transitions. The predicted values of these, according to 
Eq. 1271 are shown as horizontal line segments near the stress axis. 




FIG. 3: Model I Flow curves for a = 0.95, with U2 > f 2 = 4. Steady state flow occurs only for 
stresses above the yield stress. 



stress is increased beyond a yield stress which itself increases monotonically with V2- (No 
new features arise at larger V2 than those shown here.) The application of stresses beyond 
the yield stress lead to 'shear melting' of the glass, followed by a thinning viscosity, similar 
to the behaviour at stresses above ac2 in Fig. |21 

Our results for Model I with a = 0.95 can then be summarised as follows: far from 
the quiescent glass transition {v2 < 3), the rheology is quasi- Newtonian; for 3 < f 2 < 3.8, 
shear thinning is apparent. In the range 3.8 ;^ ^2 < 4, shear thickening becomes prominent: 
upon increasing V2 the flow curve first becomes nonmonotonic, and then (for V2 ^ 3.866) 
the curve reaches all the way back to the stress axis, indicating a jamming transition to a 
nonergodic state. Finally, for f 2 > 4, a conventional yield stress appears, below which there 
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FIG. 4: Model I flow curves at V2 = 3.9 for various a. For a = 0.9, jamming transitions and their 
analytical locations are indicated as in Fig. \^ 



is no (steady-state) flow. At stresses exceeding this yield stress, behaviour is similar to that 
above ac2 in the 'full jamming' flow curves. 

3. Dependence on a 



We now fix the value of V2 and investigate the effect of varying a. As mentioned earlier, 
the experimental meaning of a is not obvious: it controls the extent to which the mode 
coupling vertex increases with applied stress, and seems related to the susceptibility of the 
static structure to external forcing. Here we choose V2 = 3.9 (for which the zero-shear 
viscosity is r ^ 814): this is close enough to the quiescent glass transition to uncover some 
interesting behaviour. This behaviour arises for a values of order unity; note that if instead 
a <C 1 there is no shear thickening. Flow curves for a range of a < 1 are shown in Fig. |3 
The variation in the curves with increasing a is qualitatively similar to that in Fig. |21 with 
the shear thickening becoming stronger as a is increased. 

In fact there is an important qualitative distinction between Fig. E] and Fig. |21 which is 
not immediately apparent from the plots. Specifically, as a is varied, the asymptotic slope 
of the flow curves at large stresses diverges as r ~ 1/(1 — a) (Sec. IIVA2|) . In contrast, 
varying V2 at constant a < 1 leads always to curves of the same (finite) limiting slope. For 
a > 1, the slope of the flow curve at large stresses remains infinite which means there is no 
upper branch at large 7 as illustrated in Fig. El This is consistent with our earlier conclusion 
that, for a > 1, limg-^oo t = 00. However, the topology of Fig. El is not the only possibility: 
there is a small region of parameter space (1 < a < 1.0032 for V2 = 3.9) for which the 
system refluidises upon increasing the stress, followed by a second window of nonergodicity, 
which persists indefinitely upon further increasing a. This somewhat peculiar behaviour is 
illustrated in Fig. El 

A remark is in order on the physical meaning of flow curves that show no liquid-like 
branch at high stresses. We do not expect that applying an arbitrarily large stress will fail 
to make a colloidal system flow at all; but the flow could be unsteady, or could be spatially 
nonuniform involving a fracture mechanism, for example. Neither outcome is addressable 
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FIG. 5: Model I flow curves at V2 = 3.9 with a > 1. There is no refluidisation from the jammed 
state to a high stress fluid. 
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FIG. 6: Model I flow curve for a = 1.002 (in a log-log plot), showing two windows of jamming 
Fluidisation from the second window by increasing the stress is not possible. 



within the modelling framework developed here, so that flow curves such as those of Fig. UJ 
need not be regarded as unphysical in our context. On the other hand, whenever a > 1 such 
flow curves arise even for small V2, far from the quiescent glass transition. This would surely 
be unphysical so we assume that if a ever exceeds 1, it does this only at high concentrations. 
A final, intermediate flow scenario is found at the critical value a = 1. In this case it is 
difficult to resolve some important details of the flow curve numerically: there is nontrivial 
structure in regions where a and r are large. There is certainly a 'low stress' window of 
nonergodicity for 0.165 ^ cr < 5.935 (for V2 = 3.9). Upon increasing the stress further, 
the jammed state yields, followed by a second regime of thickening in which r becomes too 
large to track numerically. This suggests that a second nonergodicity transition might occur: 
however, we shall see in Sec. IV Bl (by examining Eq. ITTj) that this is not the case. 
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FIG. 7: Phase diagrams for Model I at various V2- The curves denote jamming transitions. Ah 
states below the curve for a given value of V2 are fluid states, whilst those above (and on) the line 
are nonergodic, jammed states. 
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FIG. 8: Phase diagrams for Model I at smaller V2- (Note the different scale to Fig. d) As before, 
jammed states lie on or above the line for a given value of V2- 



B. Phase Diagrams 

The above concludes our survey of flow curves for Model I, which has covered the basic 
rheological scenarios found. To get an overview of the behaviour over a wider range of 
parameters, we now concentrate on the jamming transitions and their locations in parameter 
space. By solution of Eq. EZI for a range of f 2 and a, we can calculate 'phase diagrams' 
depicting the boundary in parameter space between ergodic and nonergodic regions. 



1. The {a, a) representation 

Phase diagrams of Model I, in the {a, a)-plane, are shown in Fig. [7| and Fig. |H| This 
representation of the phase behaviour shows that, at a given value of V2, nonergodic states 
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appear provided a is sufficiently large. As V2 is lowered below 4 (so that the system at rest 
is ergodic), the minimum value of a necessary for nonergodicity under stress increases from 
zero until, for V2 less than a critical value V2 (found below) there is no transition unless 
a > 1. Once this condition applies, a transition occurs regardless of any further decrease of 
V2 (Eq. EHl). 

In Fig. m (respectively Fig. El), we showed examples of Model I flow curves in which a 
jammed state could (could not) be refiuidised by further increasing the stress. Figs. I7|8I 
show that such refiuidisation occurs within a window of a which disappears when V2 falls 
below a second critical value V2 (also found below). Outside of this window, the upper, 
refiuidised branch of the flow curve is absent (Fig. EJ. When such a window of a is present 
at a > 1 (visible for the lower two curves in Fig. |H)) an upper branch of the flow curve exists, 
but its relaxation time diverges in the large stress limit and so this branch shows a second 
jamming transition with no ultimate fiuidisation at high stress (as in Fig. ^. Flow curves 
showing full jamming with a single refiuidised branch (as visible in Figs. I!2I4|) thus require 
both a somewhat less than 1 and V2 somewhat less than 4. 



2. The (cr, V) representation 

For any given parameter set, the (a, a) representation of the phase diagram provides 
a useful picture of the behaviour under increasing the stress. However, the full range of 
jamming scenarios is most easily summarised using a different representation. This involves 
plotting the locus of transitions in the (cr, "K) plane, where we recall that V = ^2(0"). We 
refer to the resulting transition curve as Vc{<j); it is found by solving Eq. |2Z|for V2 at fixed 
a = Oc-, and has no dependence on parameters such as a. Now, in Model I, V{(y) is just 
a straight line with slope a which intercepts the V axis at V2- Each transition point, for 
given V2 and a, corresponds to an intersection of V^(o") with Vc(cr). The range of behaviour 
obtained upon varying f 2 and a simply corresponds to the different ways in which a straight 
line may intersect the fixed locus of possible transition points represented by the curve K:(c"). 
This is illustrated in Fig. [HI 



3. Implications for V2 , V2 

Above, we stated that, for V2 < V2 multiple transitions are no longer allowed. The {a, V) 
representation allows a precise calculation of V2, as follows. From Eq. |27l the maximum 
slope of K(cr) is V^{a^ax) = 27 - IdVs ^ 1.019 at ama. = 1 + 2/^3. For a > K'(^max), 
there is then a single jamming transition if V2 < 4. For refiuidisation to occur, the jam- 
ming transition must occur at a stress a < amax'- beyond this point the slope V^'(cr) is 
monotonically decreasing and so a straight line cutting the transition curve (from below) 
beyond amax cannot encounter it subsequently. The critical value V2 is then found by setting 
V2 + V^{amax)o'max = Vc{<Jmax)', that is, V2 is the intercept on the V axis of the tangent 
to the point of maximum slope. This gives V2 = 9(1 — l/v3) ~ 3.804. This geometrical 
interpretation is illustrated in the inset to Fig. El 

The slope of the transition line Vc{cr) approaches unity at large stresses ensuring that, 
for a > 1, the system is nonergodic in the limit. In contrast, for a < 1, jamming transitions 
occur only if f2 is sufficiently large. The case a = 1 is of particular interest, as we were unable 
to elucidate its behaviour numerically. The minimum value of f 2 for which a transition occurs 
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FIG. 9: Phase diagram in the (cr, y)-plane. Above and on the transition boundary, states are 
nonergodic. The dotted Une corresponds to V2 = 3.9, a = 0.95: the two intersections of the transi- 
tion boundary correspond to jamming, followed by refluidisation. The dashed line (corresponding 
to V2 = 3.8, a = 1.1) has a single intersection, showing that for these parameters there is no 
refluidisation of the jammed state. In order to bring out the structure present in Vcicr), in the inset 
we plot Vc{cr) — a (solid line) as a function of a. We show the tangent to this curve at its point of 
maximum slope (dashed line), whose intercept with the vertical axis identifies V2- (See main text: 
the geometrical interpretation of V2 holds equally well in this representation as the subtraction of 
a linear term does not change curvature.) For a = 1 the line V{a) — ex is independent of a and has 
the value V2- Thus we identify V2 (see text) as the smallest value taken by Vc{a) — a: the dotted 
line corresponds to V — a = V2- 



(with a = 1) is 'U2, as defined earlier. This is now obtained by asking what is the minimuni 
intercept (with the \^-axis) for which a straight line of unit slope cuts the transition line 
Vc{a). This yields a value of 'U2 = 1 + 2^/2 ^ 3.828. (Again, the geometrical interpretation, 
suitably altered, is illustrated in the inset to Fig. M) At a = 1, there is no transition below 
this value of V2', for V2 < V2 < 4: there are two transitions (jamming and refluidisation at Ud 
and 0"c2 respectively). As f2 — > (1 + 2-\/2)"'", ad and ac2 merge, and as f2 — > 4~, dd -^ 
and (Tc2 — > 00. 

We have shown that our critical values of V2 obey V2 > V2- For V2 between these two 
values, flow curves with two separate jammed states (Fig. ^ can be found, but there are no 
flow curves with a single jammed state refiuidised at high stress, as seen in Fig. |21or|l| 

4- Yielding of a quiescent glass 



Finally, this representation of the phase diagram allows us to elucidate the yielding of 
a quiescent glass within Model I. As a — ;> 0, Vc{a) — > 4, and V^{a) — >• 0. Therefore, for 
a = 0, the yield stress is zero at the quiescent transition (^2 = 4) and increases smoothly 
u pon moving deeper into the g lass. This differs from the microscopic MCT calculations 



u pon moving aeeper into tne g lass, inis aiiiers irom tne microscopic MUi calculations 
of iFuchs and Gated ( 2002L 120031 ) where, although there is no explicit stress dependence of 
the vertex, the yield stress rises discontinuously from zero at the static glass transition. 
However, in Model I any a obeying < a < 1 restores this expected discontinuity of yield 
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FIG. 10: Model I correlation functions for stresses approaching ad = 0.2 for V2 = 3.9, a = 0.95. 
Sufficiently close to the transition, there is a clear separation between the /3-decay and the flow 
induced final decay, resulting in a well defined plateau. At lower stresses, the relaxation time r 
become comparable to that of the /3-relaxation, and the plateau is lost. 



stress at the static F2 transition. The resulting yield stress increases with a and diverges for 
a > 1, beyond which, within Model I, the glass cannot be shear melted at any finite stress. 



C. Relaxational dynamics 



The results presented thus far are determined solely by the relaxation time (equivalently, 
viscosity) for given parameter sets. This is the time integral of the correlator. We now 
consider the functional form, rather than just the time integral, of the correlator (f){t) in 
Model I. 

Fig. ITUl shows the correlation functions on approaching a typical jamming transition. At 
early times the correlators for a range of applied stress are very similar: the /3-relaxation 
timescale (on which the initial decay of correlations occurs) is approximately constant. The 
correlators differ significantly only in their terminal relaxation times r ~ I/7. This shear- 
induced relaxation time diverges as the jamming transition is approached (7 — > 0). The 
consequences of this are as discussed in Section IIV Al On timescales t <^ I/7, Model I 
behaves as a static F2 model with coupling constant V = V2 + ao which, unless V = 4, has 
a finite rg. Fig. ^2 illustrates the equivalence of correlations in Model I, at times t -C I/7, 
with the 'underlying', noncritical, F2 model. 

In calculating the whereabouts of transition points (Sec. IIVA|) we assumed that, in the 
vicinity of a jamming transition, correlation functions were well described by an exponential 
at late times, and cited numerical evidence. Fig. ^2 demonstrates this for one particular 
parameter set. A similar agreement has been found near jamming transitions for other 
parameters, including transitions at both o"ci and 0"c2. 
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FIG. 11: The correlation functions for Model I with V2 = 3.9, a = 0.95 and a = 0.199, compared 
with the correlation function for the F2 model with V2 = 3.9+ (0.95 x 0.199) = 4.08905, and (at late 
times) with the exponential consistent with Eq. 1201 At o" = 0.199, the shear rate is 7 = 1.04 x 10~^ 
(the initial jamming transition for these parameters is at o"ci = 0.2). 



VI. MODELS II V 

Having studied Model I in some detail, we turn now to consider results from Models 
II-V, focusing on qualitative differences between these and Model I. By doing so, we hope 
to clarify which results depend on particular choices made in setting up that model, and 
which are more robust. 



Model II 



The analytic work of Sec. II V Al applies to Model II, and so the phase behaviour is 
again determined via Eq. |23 (As with Model I, numerics show that the late-time decay of 
correlations is well described by an exponential verifying the analysis leading to Eq. 1771 ) 
Using the (a, V) representation of the phase diagram we can elucidate the phase behaviour 
of this model. The transition line Vc{(t) remains unchanged from Model I, whereas V{a) 
is now quadratically increasing with stress. Since the maximum slope of V^(cr) is finite, for 
sufficiently large stresses V{a) > Vc{cr) and so, provided only that a > 0, Model II exhibits a 
jamming transition irrespective of other parameters. This ensures that, in the limit of high 
stress. Model II is nonergodic. Stress-induced jamming is never followed by refiuidisation: 
and while a quiescent glass may be shear melted (depending on a), this will always jam 
again at higher stresses. These scenarios are illustrated in Fig. ^1 We have obtained flow 
curves (not shown) corresponding to both of these scenarios. 

The flow curves of Model II also differ from Model I at low stresses for a system that is 
liquid at rest; here one finds (for suitable parameters) a region of downward curvature near 
the origin. Hence there is an initial shear thinning region before the shear thickening and 
arrest scenarios are encountered. Fig. IT^ shows the fiow curves for Model II with V2 = 3.9, 
a = 0.2, for which full jamming is present: compare this to the fiow curves for V2 = 3.9 in 
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FIG. 12: The (cr, F) representation of the phase diagram for Model II. As in Fig. |^ intersections of 
the transition line V^(o") with V{a^ mark the jamming transitions. The quadratic stress dependence 
of V^(o") in Model II ensures at least one jamming transition. At high enough stress, Model II is 
nonergodic, although (as shown) a static glass may be shear melted before jamming. 
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100 Y 

FIG. 13: Model II flow curves for v^ = 3.9, a = 0.2. Note the presence of thinning at low stresses, 
followed by thickening and jamming. No fluidisation after stress-induced jamming is found for 
these (or any) parameters with this model. 



Model I (Fig. |^, in which no thinning is apparent prior to thickening. 



B. Model III 



We have just seen that in Model II the system always jams at high stresses but is shear 
thinning at low ones. In Sec. II I II we argued that, around a = 0, the quadratic form (Model 
II) has a sounder basis (^2(0") should be a symmetric function) whereas there is no reason 
to restrict attention to a quadratic at larger stresses. The wider variety of flow scenarios 
emerging from Model I at high stressses suggests some utility for a model in which a quadratic 
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FIG. 14: The Model III flow curve with V2 = 3.9, a = 0.95. Shear thinning preceeds thickening 
and jamming. 




FIG. 15: The Model III phase diagram in the {a, V) representation for V2 = 4.1, a = 0.8. A shear 
melted glass (or quiescent fluid; not shown) can exhibit thickening and full jamming at larger stress, 
as in Model II (Fig. I12|) . But in this case, the ultimate state at high stress is a fluid for a < 1. 



stress dependence at low stresses goes over to a less severe form at higher ones. Model III 
(see Eq. E} provides a simple way to allow for this. 

As shown in Fig. ^^ in Model III shear thinning can preceed thickening, just as in Model 
II. Unlike Model II however, full jamming may now be followed by refiuidisation: the high 
stress limit is ergodic if a is smaller than unity. This is illustrated by the phase diagram 
trajectory shown in Fig. UHl in which a quiescent glass is first fiuidised, then jams, and is 
finally refiuidised again. 



Model IV 



M odel IV incorporates a n onexponential form of f{'jt) inspired by the full MCT calcula- 



tions 



Fuchs and Gated (J2002)]. This has not been treated analytically, and so results were 



obtained solely by numerical means. We find that the rheology of this model is qualitatively 



24 



Numerically calculated correlator 
Exponential fit to tail 



Log(t) 



FIG. 16: The Model IV correlation function for V2 = 3.9, a = 0.3 and a = 0.78 (the system has 
been refluidised: ac2 ~ 0.775), with an exponential fit to the late-time decay. The late-time decay 
deviates from exponential form. 



similar to that of Model I: the same shear thickening scenarios appear, and so we do not 
show any flow curves. 

In contrast, the correlation functions do differ from those of Model I. In particular, the 
late-time exponential decay present in Models I and II does not consistently appear in this 
case: in Fig. ^J we show an example where the late-time correlator is not well described 
by an exponential. This suggests that the exponential decay of Models I, II is an artefact 
of their precise form (nevertheless, it is helpful as it allows analytic progress). A second 
difference regards the behaviour in the case a = 0. In Sec. |3 we showed that, for Model 
I, the a = yield stress vanishes smoothly as V2 -^ 4+. Here, elucidating this point is not 
so simple, as we have no analytic expression for the transition points. However, numerics 
suggest that the yield stress does not vanish smoothly as V2 -^ 4+ for Model IV (Fig. ITTjl . 
We have no explanation for the very small flow rates needed to reveal this effect. 

D. Model V 



We turn finally to Model V, which differs from Model I in that it is based on the F12 
schematic model, rather than the simpler F2 model. This creates an extra dimension in 
the parameter space, since there are two coupling parameters whose variation with applied 
stress need not be identical - ie, in Eq. ^1 ai need not be identical to a2- Our exploration 
of the parameter space has uncovered thickening and jamming scenarios qualitatively similar 
to that of Models I and IV. Other scenarios may be present, but we have not found any 
which differ qualitatively from those shown earlier. 

However, as in Model IV, the correlation functions do not always show a late-time expo- 
nential decay in the vicinity of jamming transitions. An unresolved issue concerns nature of 
the a = yield stress: does it appear discontinuously at the quiescent glass transition (as 
in Model IV), or does it vanish smoothly (Models I, II and HI)? We have been unable to 
answer this question numerically: but in contrast to Model IV, the yield stress, if it exists 
at all, is very small (< 10~^). 
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FIG. 17: The a = flow curve of Model IV at V2 = 4. Note the logarithmic scale of the 7 axis. At 
sufficiently small shear rates, the flow curve levels off, and is constant over fourteen or so decades 
in the shear rate, suggesting that in this case the yield stress does appear discontinuously at the 
glass transition. Whilst the shear rates required to capture this behaviour are exceedingly small, 
the yield stress itself is not. 



E. Structural stability 



Having described the results of the various model variants, we now discuss the structural 
stability of our class of models. The shear thickening/jamming scenarios that we find appear 
fairly robust with respect to the precise form of the memory function, indicating that they 
arise generically from the competition between stress induced jamming and flow induced 
fluidisation, rather than from the precise mathematical forms chosen here. The exception to 
this is the case of Model II (and, more generally, any case in which the coupling parameters 
increase too quickly at large stress) in which fluid states are precluded at high stress regard- 
less of proximity to the quiescent glass transition. This is unphysical: very dilute colloids 
do not jam under flow. However, Model II is still useful in that it takes the 'correct' form 
around o" = 0, suggesting that thinning can preceed shear thickening. Model HI successfully 
combines this with the desirable features of Model I, in which fluid states at high stress are 
allowed, but only for a; < 1. Since a (like ^2) can be concentration dependent, this allows 
an intriguing possibility, that shear melting of colloidal glasses (stress induced or otherwise) 
is precluded above a certain threshold of concentration. 

In contrast, neither the exponential decay nor the smoothly appearing yield stress at 
q; = found in Models I - HI are robust under alterations to the chosen form of the memory 
function. In conventional MCT, late-time decays in the vicinity of the glass transition are 
not exponential, but rather are well described by the stretched exponential form (e"*^*/"^-* 
with the exponent {3 < 1). This 'dynamical stretching', which is found experimentally, is 
absent in the F2 model, although it reappears in F12 and other, more general, schematic 
models. It is not surprising, then, that Model I leads to exponential decay whilst Model 
IV does not. Dynamical stretching is associated with the co-operative motion required for 
relaxation in a dense liquid. In the vicinity of a jamming transition, relaxation occurs via 
smaller scale motions which are then accentuated by the flow. Thus, although our results 
regarding this point are inconclusive, one might expect that dynamical stretching is less 
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prono unced near these transitions than near the quiescent glass transition Fuchs and Gated 



(l2003l )]. If so, the assumption of an exponential decay could be a harmless one. 



VII. DISCUSSION 
A. The jamming transitions 

Nonergodic solutions to our schematic MCT equations are always available for V > 4 (Eq. 
[T5|) . However, they are not the only possible solutions. We have seen that ergodic, flowing 
solutions can also exist within this region: if both are present, the iteration method we used 
for T (based eg on Eq. EH) leads always to the ergodic solution, and jamming transitions 
occur at points where these solutions cease to exist. So far, we have tacitly assumed that 
physical stability follows from such stability under iteration: although plausible, this is not 
obvious. We discuss this assumption below, but firstly, presuming that it holds, we compare 
our jamming transitions to the conventional static glass transition as described by standard 
MCT. 

In standard MCT, glass transitions occur at points in paramet er spac e wher e nonergodic 
solutions become available. At these so-called 'critical points' Gotzd (l98y)], there is a 



bifurcation in the long-time limit (/) of the correlation function: it jumps from its value in 
the fluid (/ = 0) to a nonzero value fc, followed by a nonanalytic variation as the coupling is 
further increased {f — fc ^ e^^^, where e is the separation from the transition in the space of 
coupling constants). At lower values of the coupling, no such nonergodic solution is possible. 
In the case of our jamming transitions, this is not the case: rather (at fixed V2, a . . .) the 
system jams at values of the stress for which the ergodic, flowing solutions vanish. These 
points are not critical points, and do not correspond to a bifurcation in /. Therefore, the 
initial jump in / at jamming transitions is followed by an analytic variation with coupling 
(ze, with the stress if other parameters are fixed). This point is illustrated by the variation 
of the nonergodicity parameter / with a, as shown in Fig. ^[ which corresponds to one of 
the 'full jamming' flow curves in Fig. |21or|3]for Model I. 

In the figure, there is a range of stress for which an arrested solution has become available 
(indicated by the existence of nonzero / in the upper part of the diagram) but is not chosen - 
rather, the system remains flowing (as indicated in the lower part by a nonzero shear rate). 
Only when the viscosity of the flowing system diverges (at cTci), signalling that flowing 
solutions are no longer possible, does the system become nonergodic. At 0"c2, a flowing 
solution again becomes available. 

The fact that jamming transitions are not critical points was used in Sec. llVl to guide 
the analysis presented there. As discussed also in Sec. IV CI the /5-relaxation time rg, and 
not just the longer a-relaxation, diverges at a critical point in MCT, whereas at a generic 
jamming transition rg remains finite. This means that, by moving close enough to a jamming 
transition, we can make the shear rate sufficiently small that jTjj -C 1 so that the underlying 
/3-relaxation is essentially unperturbed by the flow. This gives a a pre-plateau relaxation 
of the correlator on time scale Tj3, followed at time r ~ I/7 by the flow induced terminal 
relaxation. This leads to a very well defined plateau region, visible in Fig. ^[ which is quite 
different from a quiescent MCT scenario. 

This prediction is testable in principle, although the advection of density fluctuations 
makes the measurement of correlators somewhat tricky. One idea is to monitor the de- 
cay of density fluctuations with wavevector along the neutral (vorticity) direction. These 
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FIG. 18: Schematic variation of 7 and / with a (Model I, a < 1, full jamming scenario). Solid and 
dashed lines denote stable and unstable solutions under the iteration described in Sec. IIVI Arrows 
indicate flow of 7 under the iteration. Note that the ergodic solution remains stable beyond the 
critical point a* = (4 — V2)/a where a nonergodic solution first exists, ie, V2 + ao > 4. 



fluctuations are not advected by the flow, so that there is some chance of observing the 
characteristic properties of the correlators just described, for example at stresses just above 



B. Physical stability of solutions 



We have so far assumed that stability of solutions to our MCT equations under iteration 
{eg of Eq. 1^^ is equivalent to physical stability of those solutions. There are two ways 
in which this assumption can fail. One is well understood, and relates to the mechanical 
instability of flow curves of negative slope with respect to shear banding. We discuss this 
below (Sec. IV11B2|) . but first address a different question - that of selection between 
different solutions of the schematic MCT equations when more than one exists. 



1. Iterative versus physical stability 



As mentioned above, iterative stability selects ergodic solutions over nonergodic ones, 
when both exist. Although such iteration does not directly map onto temporal evolution of 
the shear rate, one could imagine a transient violation of Eq. El resulting in an infinitesimal 
shear rate arising in an otherwise nonergodic state. Assuming a correspondence between 
iteration and time evolution, this would carry the system to the flowing solution whenever 
it exists; the arrows in Fig. ITHl would be genuine trajectories. Without knowing this cor- 
respondence, however, it seems wise to allow that when nonergodic and flowing solutions 
of our models both exist mathematically, either or both could be locally stable. For the 
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FIG. 19: Full jamming scenario of Model I, showing in bold the full manifold of steady state MCT 
solutions, including nonergodic solutions that are iteratively unstable (these lie on the vertical axis, 
outside the window ad < a < ac2)- Dotted line: see text. Vertical dashed line: extreme possible 
discontinuous shear-thickening curve arising from shear-banding. Bold arrows: extreme possible 
hysteresis loop. Light arrows: possible exploration of further jammed states. 



full jamming scenario of Fig. ^Jand Figs. |2l El this leads to the manifold of steady state 
solutions depicted by the bold line in Fig. [THl 

A pessimistic view, on the other hand, is that the nonergodic solution is the only physi- 
cally stable one whenever it does exist (so that physical stability is the opposite of iterative 
stability). This assumption would destroy much of the interesting rheology we have re- 
ported; for example, the full jamming scenario seen in Fig. ^J would be replaced by a flow 
curve that consists of the lowest branch of the S-shaped flow curve for stresses less than a*, 
then jumps discontinuously back to 7 = (dotted line) and remains there for all a > a* 
(vertical bold line). This scenario cannot really explain either continuous or discontinuous 
shear thickening. Moreover, its presumption that a nonergodic state is always preferred to 
flowing alternative contradicts the usual notion that glasses are arrested because they are 
trapped. 

2. Mechanical instability and hysteresis 



Let us assume instead that the full manifold of steady-state solutions are candidates for 
physical stability. These solutions describe putatitive steady, homogeneous states; however 
it is well known t hat a n y flow curve con taining regions of negative slope is unstable to shear 
banding DhontI ( 19991 1. lOlmstedl ( 19991 )]. For shear thickening systems the bands comprise 
slabs of material oriented layerwise in planes perpendicular to the vorticity (neutral) axis; 
these slabs have a common shear rate but different stresses. (The shear thinning case, in 
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which lay e rs are oriented parall e l to t he sliding walls of a couette, may be more familiar 



vnicn lay e rs are orientea paraii e i to t ne snamg wans oi a couette, may De more lamniar 
Olmstedl (19991), ISpenlev et al\ ( 19931 )].) In a rheometrical experiment with solid walls. 



these stresses get averaged over the slabs; thus, even a stress-controlled experiment controls 
only this average. 

Accordingly the part of the solution manifold in Fig. ITIH between ad on the vertical axis 
and the point of vertical tangency to the flow curve (at stress ar) is mechanically unstable 
and will not be observed in steady state, even in a controlled-stress rheometer. The standard 
shear banding scenario then predicts that if a quiescent fluid is progressively sheared with 
increasing 7, at some point before the vertical tangent is reached the stress will jump onto 
the upper branch of the flow curve. This is an extreme limit of stability; in practice it could 
jump sooner, at a shear rate that could depend on variou s aspects of sample history. (Such 
phenomena are well studied in the shear-thinning case; see lGrand et al\ ( 19971 ).) Thus, shear 



banding explains discontinuous shear thickening as the natural consequence of an S-shaped 
flow curve. If 7 is increased further, the upper branch will be followed out to the right. 

If 7 is now slowly reduced, there can be hysteresis, in which the jump back to the 
lower branch is delayed. For an S-shaped curve that does not cross the axis (see Fig. |21for 
examples), this must occur at or before a second point of vertical tang e ncy, a n d the hystere sis 
loop is of limited width (as reported experimentally in Frith et al\ (1996|), iLaunl (199^ 



However, in the 'full jamming' curve of Fig. ^Jit becomes possible in principle to track right 
back to the vertical axis, achieving the jammed state of zero shear rate but a finite stress. 
This offers one way to access the jamming transition at ac2 and the stress-induced glass at 
c" < o'c2- Another is to perform a stress-controlled measurement in the window between a^ 
and ac2', however, since the controlled stress is a spatial average, this might lead to a banded 
flow instead (corresponding, eg, to a point on the vertical dashed line in Fig. ITUj) . Once the 
vertical axis has been reached, it should be possible to vary the static stress, at least within 
the window of full jamming. It should also be possible to see whether the nonergodic states 
outside this window are locally stable; if they are, the stress could be increased indeflnitely 
(in principle!) without steady flow arising, and/or reduced as far as a*. The latter sets the 
lower left corner of an extreme possible hysteresis loop, as indicated in Fig. ^| Note that 
we have entirely excluded creep from these discussions which address only steady-state flow 
and its absence; in practice this could complicate matters considerably. 

Note also that if a model is chosen with no upper branch to the flow curve, as in Fig. 
El (Model I with a > 1) the preceding discussion of hysteresis might still apply so long as 
such a branch is furnished by some other mechanism, such as fracture. Experi mentally, it 



often reported than flow in the strong shear thickening regime is far from steady [Frith et al. 



( 19961 )]: this might indeed be consistent with a flow curve such as Fig. [SJ plus an unsteady 



fracture mechanism at high stresses. This hybrid mechanism would allow a complete sepa- 
ration of stress scales between the upper (shear-thickened) part of the effective flow curve 
and the lower (quasi- Newtonian) part; in particular, the upper stress level could vastly ex- 
ceed the intrinsic modulus scale for a colloidal glass. The latter is of order k^T ja?' close 
to the arrest transition, which is of order 1-100 Pa for typical parameters and decreasing 
with particle size. M u ch higher st r esses are indeed sometimes reported in jammed colloids 
d'Haene et ^ (jl993[ ). lFrith et d\ (jl99f " 



30 



C. Comparison with experimental jamming transitions 

As discussed in the Introduction, phenome na reminiscent of our i amming transitions have 
been observed experimentally. In the work of JBertrand et ali ( 20021 ). jamming occurred only 



for a particular range of concentrations: at lower volume fractions, shear thickening was 
observed, whilst at higher concentrations, the sample was solid in its quiescent state. A 
similar sequence arises in our Model I if V2 increases with concentration (as it surely does) 
and a remains less than unity. Thus it is tempting to identify the iar nming transition s of ou r 
model with the stress-induced solids recently found experimentally Bertrand et all (|2002[ ). 



Haw (200 



However, su c h an identification is uncertain. The particles used in the experiments of 
Bertrand et all ( 20021 1 are perhaps too large for glass transition concepts to be relevant; 



they have radius a ~ 3/im, so that Brownian forces are small compared with forces such 
as gravity (no density matching was done) and a purely mechanical descripti on m i ght be 
more appropriate. The particles which most easily jam in the experiments of iHawl (J200^ 



are also at the upper end of the colloidal size range, with radii a ~ l/im; it would be of 
great interest to study shear-thickening with particle sizes smaller than this. Also, most 
experiments on dense suspensi on rh e ology d o not utilise the cl ean model systems latterly 
favoured by colloid physicists Haw! (|2004 ). IPham et al\ ( 20021 )]. It is possible to closely 



approximate hard-sphere interactions (this generally requires optical index matching with 
the solvent) to which attractive forces can then be added controUably. Further investigations 
of shear thickening in such model systems would be particularly welcome. Ideally one would 
like to know the behaviour as a function of concentration, particle size, and interactions. 
One could then aim to establish a rational correlation between these and parameters such 
as V2 and a in our models. 

VIII. CONCLUSIONS 

The rheology of concentrated suspensions is sufficiently complex that a full microscopic 
theory is not yet feasible. Deliberately simplified models can prove useful in improving our 
understanding of these materials. In this article we have presented a family of models of the 
rheological behaviour of very dense suspensions, motivated by the Mode-Coupling Theory 
of the glass transition, but neglecting hydrodynamics. 

Upon variation of model parameters, a wide range of rheological behaviour was found, 
ranging from shear thinning, through continuous and discontinuous thickening (the latter 
interpreted as shear-banding), to a jamming transition into a nonergodic solid state. De- 
pending on model parameters, this may or may not be refiuidised again at higher stresses; 
in cases where the model predicts no refiuidisation one might expect another mechanism, 
such as fracture and/or unsteady fiow, to intervene. (There is some experimental evidence 
for this at very high concentrations.) Various hysteresis scenarios involving jammed states 
and/or shear banding are also possible. 

The main parameters in our simplest model variant (Model I) are V2, which is a measure 
of the distance from the glass transition in the quiescent state (this lies at V2 =4), and a 
coupling constant a which controls, in effect, the stress dependence of this distance. There 
is also some freedom in how we model the role of strain in eroding memory, but these details 
are less important to the results. 

The jamming transitions within our models differ from static glass transitions, in that 
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the transition points are not critical points of our schematic MCT: the /3-relaxation time 
does not diverge at a jamming transition. The resulting variation of the non-ergodicity 
parameter / at jamming has a discontinuity but is analytic beyond that, unlike the MCT 
glass transition which has a square-root contribution. Finally, in Models I-III, the late time 
decay of correlators close to a jamming transition is exponential. However, as discussed in 
Sec. IVIE| this is not robust but depends on the particular form of these models. Regardless 
of this, in all the scenarios suggested in this work, we expect that this final relaxation should 
differ from that close to a conventional glass transition, since the two relaxations occur in 
quite different ways. 

Although the models presented for shear thickening and jamming remain both schematic 
and prov isional, they may help guide future progress with more formal MCT-based devel- 
opments Fuchs and Gated fj2002, .2003. )]. We also hope that this work will promote a more 



careful experimental examination of the flow behaviour of very concentrated suspensions, 
preferably in systems of small, density-matched particles (so that Brownian motion is strong 
and gravity weak) with well-controlled colloidal interactions. The possible interplay between 
the jamming we describe and capillary forces at the surface of a droplet, or at the fracture 
surface of a bulk sample ruptured by shear, remains an open topic for future study. In the 
longer term, it would be very useful to develop an MCT-like theory, even at the schematic 
level, that can address time-dependent phenomena as well as the steady-state properties 
addressed in this work. 
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